This script shows results of our analysis of NDWI trend in European peatlands.

Pie charts of positive & negative counts of trend points for each monthly (May-October)

### ------------------------------------------------------------------------------------------------------------ ###
# pie charts
### ------------------------------------------------------------------------------------------------------------ ###

#categorize trend points in stron/weak negative/positive trends

point_dens_pos=lapply(slope_filtered_sig,function(m){
  m_df=as.data.frame(m[,2:ncol(m)])
  m_df_pos3=m_df%>% filter(Sslope > 0.003)
  m_df_pos3_class=cbind(m_df_pos3, sign=rep('pos3', nrow(m_df_pos3)), mon=rep(str_to_title(m[1,1]$mon), nrow(m_df_pos3)))
  m_df_pos1=m_df%>% filter(Sslope > 0.001) %>% filter(Sslope <= 0.003)
  m_df_pos1_class=cbind(m_df_pos1, sign=rep('pos1', nrow(m_df_pos1)), mon=rep(str_to_title(m[1,1]$mon), nrow(m_df_pos1)))
  m_df_pos=m_df%>% filter(Sslope > 0) %>% filter(Sslope <= 0.001)
  m_df_pos_class=cbind(m_df_pos, sign=rep('pos', nrow(m_df_pos)), mon=rep(str_to_title(m[1,1]$mon), nrow(m_df_pos)))
  
  m_df_neg=m_df%>% filter(Sslope < 0) %>% filter(Sslope >= -0.001)
  m_df_neg_class=cbind(m_df_neg, sign=rep('neg', nrow(m_df_neg)), mon=rep(str_to_title(m[1,1]$mon), nrow(m_df_neg)))
  m_df_neg1=m_df%>% filter(Sslope < -0.001) %>% filter(Sslope >= -0.003)
  m_df_neg1_class=cbind(m_df_neg1, sign=rep('neg1', nrow(m_df_neg1)), mon=rep(str_to_title(m[1,1]$mon), nrow(m_df_neg1)))
  m_df_neg3=m_df%>% filter(Sslope < -0.003)
  m_df_neg3_class=cbind(m_df_neg3, sign=rep('neg3', nrow(m_df_neg3)), mon=rep(str_to_title(m[1,1]$mon), nrow(m_df_neg3)))
  m_df_class=rbind(m_df_pos3_class, m_df_pos1_class, m_df_pos_class, m_df_neg_class, m_df_neg1_class, m_df_neg3_class)
  colnames(m_df_class)=c(colnames(m_df)[1:c(ncol(m_df)-1)],'geometry', 'sign', 'mon')
  return(m_df_class)
})
m_df_class_all=do.call(rbind,point_dens_pos)

m_df_class_all[,7]=as.factor(m_df_class_all[,7])
m_df_class_all[,8]=factor(m_df_class_all[,8], levels=month.name[5:10])

# count number of trend points per category
freq_out=lapply(point_dens_pos, function(m){table(m$sign) })

#create pie charts
lapply(freq_out, function(m){
  pie(c(sum(m[4:6]),sum(m[1:3])), labels = c('positive', 'negative'), col=c("#209fb6","goldenrod3"))+ 
  theme(labels.text= element_text(size = 50))
})

positive negative positive negative positive negative positive negative positive negative positive negative

## [[1]]
## NULL
## 
## [[2]]
## NULL
## 
## [[3]]
## NULL
## 
## [[4]]
## NULL
## 
## [[5]]
## NULL
## 
## [[6]]
## NULL

Calculate mean trend of summer months mean

### ------------------------------------------------------------------------------------------------------------ ###
##calculate mean trend of summer months mean 
### ------------------------------------------------------------------------------------------------------------ ###

#extract common coordinates (Jun-Aug) for points with significant trends
equal_points67=slope_filtered_sig[[2]][which(st_geometry(slope_filtered_sig[[2]]) %in% st_geometry(slope_filtered_sig[[3]])),]
equal_points672=equal_points67[which(st_geometry(equal_points67) %in% st_geometry(slope_filtered_sig[[4]])),]


Sslope_jun_common=slope_filtered_sig[[2]][which(st_geometry(slope_filtered_sig[[2]]) %in% st_geometry(equal_points672)),]
Sslope_jul_common=slope_filtered_sig[[3]][which(st_geometry(slope_filtered_sig[[3]]) %in% st_geometry(equal_points672)),]
Sslope_aug_common=slope_filtered_sig[[4]][which(st_geometry(slope_filtered_sig[[4]]) %in% st_geometry(equal_points672)),]
summer_months_bind=cbind(Sslope_jul_common,Sslope_jun_common,Sslope_aug_common)

summer_months_bind_df=as.data.frame(cbind(as.data.frame(summer_months_bind[,6])[,1], as.data.frame(summer_months_bind[,12])[,1], as.data.frame(summer_months_bind[,18])))
colnames(summer_months_bind_df)=c('Jun','Jul','Aug','geometry')


summer_months_bind_df[,1]=as.numeric(summer_months_bind_df[,1])
summer_months_bind_df[,2]=as.numeric(summer_months_bind_df[,2])
summer_months_bind_df[,3]=as.numeric(summer_months_bind_df[,3])

#calculate mean
summer_months_mean = lapply(seq(1,nrow(summer_months_bind_df)), function(r){
  mean_ft=mean(c(summer_months_bind_df$Jul[r],summer_months_bind_df$Jun[r],summer_months_bind_df$Aug[r]), na.rm=T)
  return(mean_ft)
  })
summer_months_mean_bind=as.data.frame(do.call(rbind,summer_months_mean)) 
summer_months_mean_cbind=cbind(summer_months_mean_bind, summer_months_bind_df[,4])
summer_months_mean_cbind_df=as.data.frame(summer_months_mean_cbind)
summer_months_mean_cbind_sf=st_as_sf(summer_months_mean_cbind_df, geom=summer_months_mean_cbind_df$geometry)

summer_months_mean_cbind_sf=summer_months_mean_cbind_sf[,-2]

#plot equal points summer months
summer_months_mean_cbind_sf$V1=as.numeric(summer_months_mean_cbind_sf$V1)

ggplot() +
    geom_sf(data=map, color='darkgrey',fill = 'white', size=0.05) + 
    geom_sf(data =summer_months_mean_cbind_sf, aes(colour=summer_months_mean_cbind_sf$V1), size = .1) + 
    #facet_wrap(~Month)+
    coord_sf(xlim = c(-30, 33), ylim = c(47, 71)) +
    scale_colour_stepsn(colours = pal(7), limits= c(-0.006,0.006), breaks =c(-0.006,-0.003,-0.001,0,0.001,0.003,0.006))+
    labs(colour = "SensSlope") +
    myTheme

50 ° N 55 ° N 60 ° N 65 ° N 70 ° N 30 ° W 20 ° W 10 ° W 0 ° 10 ° E 20 ° E 30 ° E -0.006 -0.003 -0.001 0.000 0.001 0.003 0.006 SensSlope

Calculate mean trend of overall mean

### ------------------------------------------------------------------------------------------------------------ ###
##calculate mean trend of overall mean 
### ------------------------------------------------------------------------------------------------------------ ###

#extract common coordinates (May-Oct) for points with significant trends
equal_points6721=equal_points672[which(st_geometry(equal_points672) %in% st_geometry(slope_filtered_sig[[1]])),]
equal_points67215=equal_points6721[which(st_geometry(equal_points6721) %in% st_geometry(slope_filtered_sig[[5]])),]
equal_points672156=equal_points67215[which(st_geometry(equal_points67215) %in% st_geometry(slope_filtered_sig[[6]])),]

Sslope_jun_common=slope_filtered_sig[[2]][which(st_geometry(slope_filtered_sig[[2]]) %in% st_geometry(equal_points672156)),]
Sslope_jul_common=slope_filtered_sig[[3]][which(st_geometry(slope_filtered_sig[[3]]) %in% st_geometry(equal_points672156)),]
Sslope_aug_common=slope_filtered_sig[[4]][which(st_geometry(slope_filtered_sig[[4]]) %in% st_geometry(equal_points672156)),]
Sslope_may_common=slope_filtered_sig[[1]][which(st_geometry(slope_filtered_sig[[1]]) %in% st_geometry(equal_points672156)),]
Sslope_sep_common=slope_filtered_sig[[5]][which(st_geometry(slope_filtered_sig[[5]]) %in% st_geometry(equal_points672156)),]
Sslope_oct_common=slope_filtered_sig[[6]][which(st_geometry(slope_filtered_sig[[6]]) %in% st_geometry(equal_points672156)),]

summer_months_bind_allmon=cbind(Sslope_may_common, Sslope_jun_common,Sslope_jul_common,Sslope_aug_common, Sslope_sep_common, Sslope_oct_common)

summer_months_bind_allmon=as.data.frame(cbind(as.data.frame(summer_months_bind_allmon[,6])[,1], as.data.frame(summer_months_bind_allmon[,12])[,1], as.data.frame(summer_months_bind_allmon[,18])[,1],
                                    as.data.frame(summer_months_bind_allmon[,24])[,1],as.data.frame(summer_months_bind_allmon[,30])[,1],as.data.frame(summer_months_bind_allmon[,36])))
colnames(summer_months_bind_allmon)=c('May','Jun','Jul','Aug','Sep', 'Oct','geometry')

summer_months_bind_allmon[,1]=as.numeric(summer_months_bind_allmon[,1])
summer_months_bind_allmon[,2]=as.numeric(summer_months_bind_allmon[,2])
summer_months_bind_allmon[,3]=as.numeric(summer_months_bind_allmon[,3])
summer_months_bind_allmon[,4]=as.numeric(summer_months_bind_allmon[,4])
summer_months_bind_allmon[,5]=as.numeric(summer_months_bind_allmon[,5])
summer_months_bind_allmon[,6]=as.numeric(summer_months_bind_allmon[,6])

#calculate mean
summer_months_mean_allmon = lapply(seq(1,nrow(summer_months_bind_allmon)), function(r){
  mean_ft=mean(c(summer_months_bind_allmon$May[r],summer_months_bind_allmon$Jul[r],summer_months_bind_allmon$Jun[r],summer_months_bind_allmon$Aug[r],summer_months_bind_allmon$Sep[r],summer_months_bind_allmon$Oct[r]), na.rm=T)
  return(mean_ft)
  })
summer_months_mean_allmon_bind=as.data.frame(do.call(rbind,summer_months_mean_allmon)) 
summer_months_mean_allmon_cbind=cbind(summer_months_mean_allmon_bind, summer_months_bind_allmon[,7])
summer_months_mean_allmon_cbind_df=as.data.frame(summer_months_mean_allmon_cbind)
summer_months_mean_allmon_cbind_sf=st_as_sf(summer_months_mean_allmon_cbind_df, geom=summer_months_mean_allmon_cbind_df$geometry)

ggplot() +
    geom_sf(data=map, color='darkgrey',fill = 'white', size=0.05) + 
    geom_sf(data =summer_months_mean_allmon_cbind_sf, aes(colour=summer_months_mean_allmon_cbind_sf$V1), size = .1) + 
    #facet_wrap(~Month)+
    coord_sf(xlim = c(-30, 33), ylim = c(47, 71)) +
    scale_colour_stepsn(colours = pal(7), limits= c(-0.006,0.006), breaks =c(-0.006,-0.003,-0.001,0,0.001,0.003,0.006))+
    labs(colour = "SensSlope")+
    myTheme

50 ° N 55 ° N 60 ° N 65 ° N 70 ° N 30 ° W 20 ° W 10 ° W 0 ° 10 ° E 20 ° E 30 ° E -0.006 -0.003 -0.001 0.000 0.001 0.003 0.006 SensSlope

Plot pie charts of summer months and overall mean

### ------------------------------------------------------------------------------------------------------------ ###
##plot pie charts of summer months and overall mean 
### ------------------------------------------------------------------------------------------------------------ ###
#######
### pie charts mean summer
summer_mean_pos_count=summer_months_mean_cbind_sf %>% filter(V1>0)
summer_mean_neg_count=summer_months_mean_cbind_sf %>% filter(V1<0)

pie(c(nrow(summer_mean_pos_count),nrow(summer_mean_neg_count)), labels = c('positive', 'negative'), col=c("#209fb6","goldenrod3"))+ 
  theme(labels.text= element_text(size = 50))

positive negative

## NULL
### pie charts mean all month
allmon_mean_pos_count=summer_months_mean_allmon_cbind_sf %>% filter(V1>0)
allmon_mean_neg_count=summer_months_mean_allmon_cbind_sf %>% filter(V1<0)

pie(c(nrow(allmon_mean_pos_count),nrow(allmon_mean_neg_count)), labels = c('positive', 'negative'), col=c("#209fb6","goldenrod3"))+ 
  theme(labels.text= element_text(size = 50))

positive negative

## NULL
### ------------------------------------------------------------------------------------------------------------ ###